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We perform a general group-theoretical study of the Kondo problem in monolayer and bilayer 
graphene around the charge neutrality point. Utilizing the group representation theory, we derive 
from symmetry considerations a family of the Kondo models for all symmetric placements with 
either 3- or 6-fold rotational axis of an impurity atom in an arbitrary orbital state. We find six 
possible classes of the partially anisotropic four-channel Kondo model. As the key result, we argue 
several possibilities to realize the regime of the dominant channel-symmetric two-channel Kondo 
effect, protected by the local symmetry and specifics of the graphene band structure. Our findings 
open prospects for the observation of the rich multi-channel Kondo physics in graphene and the 
associated non-Fermi-liquid behavior. 



Introduction and main results. The Kondo effect - the 
interaction of the local spin or orbital degrees of free- 
dom with the conduction electrons - plays a fundamen- 
tal role in a wide variety of condensed matter systems, 
ranging from quantum dots to strongly correlated materi- 
als [D [2] . The Kondo effect in graphene [3] has attracted 
significant interest [3^E5] due to its peculiar electronic 
properties and potential for the realization of the multi- 
channel Kondo effects. A number of theoretical studies 
of the quantum impurity problem has been undertaken 
for monolayer graphene (MLG) [^20j. and much fewer 
for bilayer graphene (BLG) ^51]. However, the analysis of 
possible Kondo models in graphene has not yet been per- 
formed in its full generality. In particular, the feasibility 
of the multichannel Kondo effect in graphene remained 
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FIG. 1: (Color online) Lattice structure of MLG (top) and 
BLG (bottom) and considered impurity atom (lA) (green) 
placements with either 3- or 6-fold rotational symmetry. 
The phases e"'^'^, e = e'^'^^^^, of the Bloch wave-functions 
XKA,KB{r) in MLG (top right) and XkA,kb(^) in BLG (bot- 
tom right) are shown, Xk'...{^) ~ Xk...{^)- 



a debated question: both pro [H [TH] and con [20] 
arguments have been put forward. 

Motivated by these interesting questions, in this Let- 
ter, we perform a general group-theoretical study of the 
Kondo effect in MLG and Bernal-stacked BLG around 
the charge neutrality point (CNP). Following the origi- 
nal recipe of Nozieres-Blandin we utilize the group 
representation theory [26l [22] to derive the family of the 
four-channel Kondo models for all symmetric placements 
with either 3- ro 6-fold rotational axis of an impurity 
atom (lA) in an arbitrary orbital state without appealing 
to any microscopic details. We find six possible classes of 
the Kondo models: there are three cases for the structure 
of the conduction electron channels and the impurity can 
be in either a singlet or doublet orbital state. 

Most importantly, we argue that in several cases of the 
impurity placement MLG and BLG band structure allows 
for the realization of the two-channel Kondo effect, where 
the exchange couplings for the pair of equivalent channels 
i?c± belonging to one two-dimensional (2D) irreducible 
representation (IR) are dominant. In these cases, in 
the low-energy regime, the symmetric two-channel Kondo 
model 



(o)(j''^"^" + r''^''^aS)4,EM (1) 



is realized in the orbital singlet case and the two-channel 
Kondo model 



^dblt 



E 

(2) 

with partially anisotropic (J^' — Jy' ^ J.: ) orbital 
Kondo interactions is realized in the orbital doublet case 
Here, ipE^ = {il'E^+t,i^E^+i,-ipE,-t^i'E,-iY is a spinor 
in the product of the channel {Ee±) and spin (t,i) 
spaces, Tj and Tj are the unity (7 — 0) and Pauli 
(7 — X, y, z) matrices in the channel and impurity or- 
bital doublet spaces, respectively, cr — {(Jx,<Jy,crz) are 
the spin Pauli matrices of the conduction electrons, and 
S = {Sx, Sy, Sz) are the impurity spin operators. 
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basis states of 


squares of 2D IRs 


MLG (a, Cei,) 
{a,Deh = Cei; x Cs) 


El + E2 


XBi- = -^(Xka + Xk'b), XEi+ = -^(Xkb + Xk'a), 

If \ I I \ 
XB2+ = yjUK'A - Xkb), XE2- = -^{Xka - Xk'b) 


Ei,2 X £'1,2 ^ Ai + A2 + El 

-^1,2 X -£^1.2 = Ai + A2 + El 


(6,C3„) 

(0, Ush — OSu X Us) 


A1+A2 + E 

A" _l_ A" _l_ 7?" 


Xe+ = Xka, Xe- = Xk'a 
XAi = ^ixKB + Xk'b), XA2 = ^{xkb - Xk'b), 


ExE = Ai+A2 + E 

Tp' ," w TT"''" ^ ' _l_ ^1 ' _l_ 


BLG (c,C3.) 


2£ = + 


XB"+ = XifA> Xe^-^Xk'A' Xei^+^Xk'b, Xe0-=Xkb 


ExE = Ai+A2 + E 


{c,Dsd = Cs^ X Ci) 


Eg + i5u 


Xb„+ = yjUifA + Xk'b),Xe^- = ^IXif'A + XifsJ 
XEg+ = -^(XkA - XK'B),XEg- = ^(Xk'a - Xifs) 


Eg,u X -Eg,!! = Alg + A2g + Eg 


(d,C3.) 


yli + ^2 + 


XB+ = Xk'A' xe- = XifA> 

XAi = -^{xkb + Xk'b), XA2 = -^{xkb - Xk'b) 


ExE = Ai+A2 + E 



TABLE I: Considered cases of impurity atom (lA) placement (first row), labeled by the rotational axis a,b,c,d and the local 
symmetry group G (Fig. [TJ . Decomposition of the 4D space Rc of the conduction electrons in terms of the IRs of G (second 
row) and the basis functions of IRs (third row). Decomposition of the squares of 2D IRs (fourth row), used for the construction 
of the exchange Hamiltonian Hj [Eq. ([5| and Tab. |ll|. 



Both models present considerable physical interest, 
largely due to the non-Fermi-liquid behavior exhibited 
in a number of regimes, but have proven challenging to 
realize in practice. Our findings thus pose MLG and BLG 
as promising materials for the realization of the rich mul- 
tichannel Kondo physics and the associated non-Fermi- 
liquid behavior. 

Graphene lattice and local symmetry of the impurity 
atom. MLG and BLG [3] are two-dimensional carbon al- 
lotropes, shown in Fig.[T] Both have a triangular Bravais 
lattice with the primitive translation vectors ai.2. MLG 
has PQ/mmm {D\f^, #191, hexagonal system) space 
group and D^h point group; its unit cell contains two 
atoms denoted A and B. BLG has P3ml (L*!^, #164, 
trigonal system) space group and D^d point group; its 
unit cell contains four atoms A, B, A, and B, two in 
each layer. 

In the presence of the lA, the spatial symmetry is re- 
duced to a point group G, with the center of the 
impurity being the fixed point. For the realization of the 
channel-symmetric 2-channel Kondo effect, it is neces- 
sary that G has 2D IRs, which is the case for hexagonal 
lattice if G contains either 3- or 6-fold rotational axis. 
Lower-symmetry point groups have no 2D IRs and the 
symmetric two-channel Kondo model is not feasible in 
these cases. Therefore, in this paper, we consider only 
such impurity placements that G contains either 3- or 6- 
fold rotational axis (these, at the same time, may be the 
likely adsorption sites). 

For MLG, there are two such vertical axes, denoted a 
and b (Fig.[T]). For the axis a going through the center of 
the carbon hexagon, G = Gqv, if the lA is out of the MLG 
z o — z mirror plane and G ~ D^h ~ C^y x Cs (equal to 
the point group of MLG), if the impurity center is in the 
mirror plane (here, Gg is the group of mirror reflection 
only). For the axis h going through the center of the 
carbon atom, G = G^y, if the lA is out of the MLG 
plane, and G = D^h — C'sd x C^, if the lA is in the MLG 



plane (i.e., lA substitutes the carbon atom). 

For BLG, there also are two such vertical axes, denoted 
as c and d (Fig.jl]). For the axis c going through the center 
of the hexagon of one layer and a carbon atom (B) of the 
other layer, G = Gsv for any of the vertical position of 
lA on the axis. For the axis d going through two carbon 
atoms {A and B) of the layers, G = D^d = Czv x Gi 
(equal to the BLG point group, Cj is the point group of 
inversion only), if the lA is at the midpoint between the 
layers, and G = C^y for any other lA placement on the 
axis. These cases are summarized in Tab. HI 

Impurity degrees of freedom. A generalized Kondo 
model [2 [2] describes minimal coupling of the conduc- 
tion electrons to the spin and, if present, orbital degrees 
of freedom of ground state of a decoupled lA, i.e., ne- 
glecting the hybridization to the conduction states. The 
decoupled lA has a definite electron occupancy number 
N, spin S (according to the Hund's rule, S — N/2 and 
5 = 2^ + 1 — N/2 for less- and more-than-half-filled or- 
bital with angular momentum I, respectively). In the 
crystalline environment of the lattice, the orbital ground 
state of lA belongs to one of the IRs Ri of the local sym- 
metry group G (we assume spin-orbit interactions weak) . 
In all cases we consider (Tab. [l]), G has only ID or 2D 
IRs [Sni HZj . Thus, the lA ground state is either an or- 
bital singlet \Sz) or doublet \a) (E)\Sz) if it belongs to one 
of the ID {Ri = Ai) or 2D {Ri = Ei) IRs, respectively. 
Here, Ai and Ei denote arbitrary ID and 2D IRs of the 
impurity orbital state, respectively, — —5*, . . . , 5", and 
a = ± arc the quantum numbers of the orbital doublet 
Ri = Ei. We choose the basis states |a = ±) so that they 
transform as e""^ under G^y ii Ri — E, and as e""'' and 
Q-a2iip ^jjjgj. jf ^. — E2,Ei, respectively. 

Thus, among the variety of microscopic possibilities of 
different I, N, and sequences of crystal field splittings, 
group-theoretically, there are only two different classes 
of the orbital state of the lA. 

Conduction electrons in MLC and BLG. The elec- 
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tronic band structure in both MLG and BLG [3 [Ml 
around zero doping is governed by four Bloch states 
X^(r) (Fig. [l]), M = KA,KB,K'A,K'B in MLG and 
/i = KA,KB,K' A.K'B in BLG, at two high-symmetry 



points K 



3|aiP 



(ai — and K' = — K in the BriUouin 



zone, referred to as valleys K and K' ^ with energy e = 
exactly at the charge neutrality point (CNP). The states 
in two valleys are related by the time reversal symme- 
try: X-Rr'...(i") = Xif. ..('')■ I'^ MLG, the two states per 
valley reside on either A oy B sublattice and are labeled 
accordingly. In BLG, the two states per valley reside on 
either Aot: B sublattice, located in different layers, while 
their weight on A and B sublattices vanishes. These four 
states X/i (r) form 4D IRs of the respective space groups 
of MLG and BLG. 

The vicinity of the CNP can then be described utilizing 
the k ■ p-theory expansion of the electron field operator 



MLG, BLG 



(r) = ^X;.(r)V;^<,(r) 



(3) 



in terms of the exact Bloch wave- functions at e = and 
the operators V';icr(r) that vary over large spatial scales 
(cr =t: i is the spin projection). In one valley, the electron 
motion in MLG and BLG is described by the chiral linear 
(massless Dirac equation) and quadratic Hamiltonians in 
the sublattice space, respectively, 




iBLG 



AB 



2m\p^_ 



AB 



Here, p± ^ Px ± ipy: Px,y = -^dx,y, v w 10*cm/s, and 
m ~ me/20. The Hamiltonian in the other valley is the 
time-reversal counterpart, h^, — h^*- These two-band 
models are valid at energies e <C i bellow the nearest 
neighbor hopping amplitude t w 3eV in MLG and e <C t_L 
below the interlayer hopping amplitude tj_ ss O.SeV in 
BLG. 

Classifying the conduction electrons. To proceed with 
the derivation of the Kondo models, we classify the con- 
duction electron states according to the local symmetry 
group G. It is sufficient to consider the states exactly 
at the CNP, treating -;/;pcr(r) as constants, in which case 
^MLG,BLG(-j.-) representation i?o of G. Ap- 

plying the symmetry operations of G to Xi-i{^) each 
considered case, we obtain the decomposition of Rc into 
the IRs of G and their basis functions (second and third 
rows of Tab. |l]). (In doing so, we keep in mind that the 
Pz atomic orbitals of Xfi (r) have odd parity under mirror 
reflection in the graphene plane and inversion.) Similarly 
to the impurity states |a = ±), we choose the basis states 
of the 2D IRs so that they transform as xe± ~ e*"^ un- 
der C^v and as XB2± ~ ^'''"^ and XEi± ~ e^^"^ under 
Cgi,. This convention eventually leads to the most natu- 
ral form of the orbital Kondo interactions (Tab. [TTl). The 
electron field operators can be rewritten in the new bases 



as \[/(r)MLG,BLG ^ X],,X»?(i")^»7(T(r), where rj label the 
basis states of the IRs, e.g., ry = Ai, A2, E+, E— in the 
MLG (5, C31,) case. For further use, we arrange the oper- 
ators ipr^ into the eight-component spinors ip presented 
in Tab. nil where -0^ = {ipri^-jiprjiY ■ For the cases (a, Cgu) 
and (6, Cst,) in MLG, our results for the classification of 
ipri agree with the earlier results of Refs. O [T31 US] ■ 

We find that the seven considered cases of the impu- 
rity placement fall into 3 different classes, presented as 
columns of Tab.[llj the conduction states split into (I) 
two different 2D IRs; (II) one 2D and two ID IRs; (III) 
two 2D IRs of the same type. 

Four conduction channels. At the length scale of the 
slowly varying operators "0^0- (r), the atomic impurity 
may be considered as a point object, whose degrees of 
freedom couple only to the conduction states with non- 
vanishing weight V'm(t(0) at its position r = 0. As the 
solution of the eigenvalue problem for h^^^,'^^'^ in po- 
lar coordinates shows, see Refs. [H HSl [T^ and Supple- 
mentary Material, for both MLG and BLG there is ex- 
actly one such radial channel for each of the four compo- 
nents fi of t/jfj^a-ir). By the standard "unfolding" proce- 
dure |2j, where the outgoing and incoming radial waves 
are mapped to the plane waves in s > and s < regions 
of the effective ID axis, respectively, these four channels 
can be represented as chiral ID channels with the kinetic 
energy 



Ho 



dsV^t(s)(-i5,-ef)7/;(s). 



(4) 



We will preserve the same notation for the ID fields ?A(s) 
as for the 2D ip{r) ones. Thus, for an atomic-size impu- 
rity, a 4-channel Kondo model is generically realized in 
MLG H [HJ Hg and BLG around the CNP. 

Kondo models from the symmetry approach. We now 
derive the Kondo models. As first outlined by Nozieres 
and Blandin T] and later implemented in a number of 
works for various Kondo systems 2 , 3n| , the most general 
possible form of the exchange interaction Hamiltonian 
Hj can be efficiently obtained based on the symmetry 
grounds without appealing to any microscopic model: Hj 
must remain invariant under all symmetry operations of 
the system. The most general form invariant under spin 
rotations reads 



Hj ^i;^0){JP + .r(TS)ij{0). 



(5) 



In Eq. ([5|, J''''^ are the operators in the orbital sector 
Rc<E} Ri'. they are 4x4 matrices in the orbital space Rc 
of the conduction electrons and, in the impurity doublet 
case, they also contain the orbital "isospin" operators 
T-y, J = 0,x,y, z (unity and Pauli matrices) acting in the 
space Ri — Ei of the states \a = ±). 

The operators J'''" must remain invariant under the 
orbital symmetry group G |31j . Such invariant form is 
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"■J 


(I) MLG: (a,C6„/D6h); BLG: (c, Dgd) 


(II) MLG: (6,C3„/i33h); BLG: (d.Cg.) 


(III) BLG: (c,C3„) 


sglt.: 
dblt.: 


\ U J 1 

( jSclSol Q \ 


h.C 

V h.C 


^ J-^i-^i (0,0) \ 
J'*^'*^ (0,0) 
(0,0)' (0,0)* J^^To / 


/ jE^'E" jE^E^ \ 

V J J J 

/ jE^E" jE^eI^ \ 
51^7 ( jB^b". jB^B-S I^^tT't 




= (V'Ecl + >V'i5„i-,^/'B,2+.'0Ec2-)' 




= (l^B° + ,V'B°-,V'B/3+,'l/'BfJ_)* 



TABLE II: The most general forms allowed by symmetry of the exchange interaction Hamiltonian Hj [Eq. ^] in MLG and 
BLG around the CNP for all possible symmetric placements with either 3- or 6-fold rotational axis of the impurity atom 
(Fig. [l] and Tab.[l]|. The entries are the expressions for J''''^; both have identical structure in but their own sets of coupling 
constants J''/ ' and J°^; ", we suppress the indices p,a for brevity. There are three classes, (I), (II), and (III) (columns), of 
the conduction electron orbital states Rc and two classes, singlet (sglt., Ri = Ai) and doublet (dblt., Ri — Ei) (rows), of the 
impurity orbital state Ri. In (I), i5ci,o2 = £1,2, £"2, £"9,14; in (II), T± = Tx ± iTy. In the orbital doublet case, the summation 
goes over 7 = 0, a;, j/, z and everywhere J^' ~ Jy" , since TxTx + TyTy is an invariant. 



efficiently constructed using the algebra of the group rep- 
resentation theory [211 [17] as fohows P'','^. The operators 

JP'" - (i?c X rI) X (Ri X Rj) (6) 

transform as a product of four representations of G, 
where Rc x i?J and Ri x Rj describe the transformation 
properties in the conduction electron and impurity sub- 
spaces, respectively. The decomposition of the product 
^ into IRs can readily be calculated ^221 [27] . The only 
allowed terms in JP'"^ are the invariants, which transform 
according to the unity IR Ai/A'i/Aig of G; each invariant 
may enter JP''^ with its own coupling constant J. The 
explicit form of these invariants is constructed by utiliz- 
ing the transformation properties of ip^a- and \a = ±) 
under G. 

This procedure yields the most general forms allowed 
by symmetry of the Kondo exchange interaction Hamil- 
tonians Hj [Eq. presented in Tab. details of the 
derivation are provided in the Supplementary Material. 
The contents of Tab. |ll) along with Tab. |l] and Eqs. ^ 
and ([5]), constitute the central result of our work. They 
describe the family of the four-channel Kondo models 

H^Ho + Hj 

in MLG and BLG in the vicinity of the CNP for aU pos- 
sible seven cases (Tab. |l]) of the symmetric placements of 
the lA with either 3 or 6-fold rotational axis. In the rest 
of the paper, we discuss their key properties and physical 
implications. 

Main properties. We find 6 possible classes of the 
Kondo models: there are 3 above-mentioned classes (I), 
(II), and (III) (columns of Tab. |ll]) for the conduction 
electron states and two classes of the impurity orbital 
state, singlet (sglt., Ri = Ai) and doublet (dblt., Ri — Ei) 
(rows of Tab. |ll|. All cases of the impurity placement 
within one class have identical structure of the Kondo 
model. 



In the orbital-singlet class (I) and (II) models, when 
the conduction states R^ break into IRs of different type, 
the channels are not coupled by the exchange interac- 
tion and the two channels belonging to the same 2D IR 
are characterized by the same coupling constant, as pro- 
tected by symmetry. These are prerequisites for the re- 
alization of the two-channel Kondo effect. On the other 
hand, in the orbital-singlet class (III) model, the conduc- 
tion sea Re = 2E = E^+El^ consists of two 2D IRs of the 
same type E (the labels a,/3 are used to distinguish be- 
tween the two subspaces). As a result the "conversion" 
processes ^ E^ , whereby the conduction electrons 
are transferred between two IRs, are present. 

In the orbital doublet case, for any class (I), (II), or 
(III) , the two channels belonging to each 2D IR Ec couple 
to the impurity via anisotropic orbital Kondo interaction 
Jq TqTq + J-^-T^Tz + J^'{txT^ + TyTy), where are the 
unity (7 — 0) and Pauli (7 = x,y,z) matrices acting in 
the space of the -Eo± states. The class (I) model consist 
of two decoupled two-channel contributions of this kind 
and there are no conversion processes Ed -H- Ec2 between 
them. In class (II) model the conversion processes Ai o 
A2 and Ai^2 ^ E between all IRs are present. In class 
(III), the conversion processes E" o E^ between two 
2D IRs of the same type are present as well, with the 
structure of the orbital Kondo interactions above. 

We note that since the kinetic energy (|4| possesses 
SU(4) rotational symmetry in the channel space i?o, the 
orbital-singlet class (III) model can be transformed to 
class (I) model without conversion processes E"' o E>^ 
by diagonalizing the matrix J'"^ (v, u' = E°',Ef^) of the 
coupling constants. Such transformation is, however, not 
possible in the orbital doublet case, since this would re- 
quire diagonalizing 3 matrices Jq'^ , J^"^ , J'[^ , at the 
same time. So the conversion processes E" o E^ can- 
not be eliminated in the orbital-doublet class (III) model. 

Increase of symmetry. We point out the following 
interesting properties. In MLG, the local symmetry 
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is increased as (a, Ce,,) {a^D^h = C^y x Cs) and 
(b^Csv) — )■ {b,D3fi = Cat, x C^) when the lA is moved 
from the out-of-plane to in-plane position. This, how- 
ever, does not lead to qualitative changes in the Kondo 
models. In the former case, the channel space i?c = 
Ei+E2^ E'{ + E'^ still breaks up into two different 2D 
IRs and still the class (I) Kondo model is realized, given 
by the sum of two 2-channel contributions with generally 
different couplings. This implies, in particular, that the 
claims [H [TH] of the 4-channel symmetric Kondo effect in 
the (a, D^h) case are unjustified: group theoretically, this 
would require the presence of a 4D IR, which is absent 
in the D^h group). In the latter case, the channel space 
= A1+A2+E A'l+A'^^E" still breaks up into two 
ID and one 2D IRs and still the class (II) Kondo model 
is realized. On the other hand, the increase in symmetry 
(c, Czv) (c, D^d) upon placing the I A halfway between 
the BLG sheets does lead to a qualitative change: the 
channel space = 2E Eg + E^ breaks up into two 
different IRs instead of two IRs of the same type] as a 
result, conversion processes vanish and the Kondo model 
is transformed from class (III) to (I). 

To summarize this part, our group-theoretical analy- 
sis allow us to make definitive symmetry-based model- 
independent conclusions about the structure of the 
Kondo model, without relying on any specific microscopic 
details. In particular, it tells exactly if and how different 
conduction channels are coupled by the exchange inter- 
action. Tab. [nj We emphasize that generally the channel 
states x,,(r) can be either pure valley states or mixtures 
thereof, as seen from Tab. |I] This analysis also resolves 
the concern [20| that "valley mixing" by the impurity po- 
tential could be detrimental for the multichannel Kondo 
effect: while valleys can indeed be mixed, the local sym- 
metry dictates that properly hybridized valley states act 
as independent channels in several instances (Tab. [lT| , 
and are also completely equivalent if they belong to one 
2D IR. 

Feasibility of the channel- symmetric two- channel 
Kondo effect. We now discuss the implications of our 
results for the realization of the multichannel Kondo ef- 
fects. It is well established [TJ [5] that the low-energy be- 
havior of the multi-channel Kondo model is determined 
by the channel(s) with the largest exchange coupling. 
Thus, effectively, the regime of the multi-channel Kondo 
effect can be realized only if the common coupling of sev- 
eral channels (belonging to the same multi-dimensional 
IR) exceeds all single-channel couplings. This condition 
proved to be extremely hard to achieve in practice: in a 
"typical" band structure with no special properties, some 
single-channel coupling will usually prevail. 

Our findings suggest that, owing to the peculiarities 
of graphene band structure, the regime of the dominant 
two-channel Kondo effect is feasible in MLG and BLG in 
several cases. 

The class (I) Kondo model, which includes the cases 



{a,Cev/DQh) in MLG and (0,^3^) in BLG, consists of 
two decoupled symmetric two-channel contributions with 
their own couplings. Whichever couplings are greater, 
the low-energy behavior will be dominated by the two 
channels of that 2D IR. 

In the BLG (d, C^v) case, if the lA is placed above the 
center of the carbon hexagon of one layer, as shown in 
Fig. [1] the atomic orbitals of the 2D IR E channels, lo- 
cated on the A sublattice in that nearby layer, are much 
close to the impurity than those of the ID IR A12 chan- 
nels, located on B sublattice in the remote layer. Thus 
the hybridization with the E channels is likely to be con- 
siderably greater than with Ai^2 channels. This should 
result in a greater exchange couplings J^^ > J-^i. 2-41,2 
and the low-energy behavior will be dominated by the 
two E channels. 

In these cases, in the low-energy regime, the symmetric 
two-channel Kondo models ([T]) and ([2| will be realized in 
the orbital singlet and doublet cases, respectively, in the 
subspace of the two dominant channels Ec±. 

Conclusion and Outlook. In summary, motivated by 
the prospect of realizing multichannel Kondo effects, we 
performed a general group-theoretical classification of the 
Kondo models in MLG and BLG in the vicinity of the 
CNP for all placements of the lA with either 3- or 6-fold 
rotational symmetry. We found six possible classes of 
the four-channel Kondo models, summarized in Tab. |IT] 
We argued several possibilities for realizing the channel- 
symmetric two-channel Kondo effect, described by the 
models ([1]) and ([2|, which are known to exhibit non- 
Fermi-liquid behavior in a number of regimes but have 
proven challenging to realize in practice. Our findings 
thus open prospects for the observation of the multi- 
channel Kondo physics in MLG and BLG, which could 
be pursued experimentally using the local probes, such 
as scanning tunneling microscopy, or transport measure- 
ments. 
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SUPPLEMENTARY MATERIAL: RADIAL CHANNELS IN GRAPHENE 



In MLG, the solution of the eigenvalue problem h^^'~^ilj — eip, X — K , K' , in the polar coordinates r = r(cos (p, sin ip) 
is provided by the eigenfunctions 

AKjAr) = M f . 'I'' ] , AK'.Ar) = r.KjAr) (MLG) (7) 

y isgnesgnj;, J|j-^+i|(fcr)e'^^^^2;'^ J 

Here, J..(fcr) are the Bessel functions and k — \e\/v. The eigenstates are characterized by the energy e, valley 
X — K, K' , and angular momentum jz quantum numbers; the latter takes half-integer values, jz = ±|, . . .. 

In BLG, the solution of the eigenvalue problem h^'^ili = 6-0 in the polar coordinates is provided by the eigenfunc- 
tions 



\^ -sgneJ|,^_i|(fcr)e'(J- /as 
where jej = ^ and the angular momentum jz — 0, ±1, ±2, . . . takes integer values. 



(8) 



The normalization of the states as 

'" ,,t 



d^rV'I,,^(r)^,,v,;(r) = 2^5(6 - e')5xx'5,.r. 



yields for the constant 



{ kl (MLG) 
'^^^'^ is the density of states, 

t (BLG) 

The expansion of the electron field operator in terms of the annihilation operators c^xj^a of the states ([t]) or ([s]) 
reads 



/I 
^ Axj^ (r)c,A,.. {X^K,K', a =t, i), (9) 
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The operators tp\a{Y) are spinors in the sublattice space, AB for MLG and AB for BLG, which are convenient to join 
into the eight component spinor 



V)(r) = 



in the valley-sublattice-spin space. 
The kinetic energy has the form 



The most general form of the exchange Hamiltonian, invariant under spin rotations, reads 



(10) 



(11) 



where JT''''^ are 4x4 matrices in the valley-sublattice space (containing the impurity orbital isospin operators in the 
doublet case). Only the states ([t]), ([s]) with nonvanishing weight at r = couple to the point impurity. This is the 
case only for jz = ±5 states in MLG, 



and only for = ±1 states in BLG, 



AB 





i sgn esgn A 



(MLG), 



AB 



V'.Aj.=-i(o) = A/; 







AB 





-sgne 



(BLG). 



AB 



Importantly, for given j^, the wavefunctions of these states reside on either one of the sublattices at r = 0. Thus, 
in both MLG and BLG, there is exactly one radial channel per valley-sublattice, and one can identify the angular 
momenta and sublattice degrees of freedom: 



3, = ]^^ A and = o B (MLG), 



(12) 



iz = -1 O i and = 1 o J5 (BLG). 
Substituting the expansion ^ into Eq. ([11]), we obtain 



where 



J'''^(e,e')=^0^^'''" 

and we join the operators of the eight states (for a given e) that couple to the impurity into the spinor 



/ 



CeKAa 



\ 



I 



CeKBa 
CeK'Aa 
\ CeK'Bcr ) 



in MLG and 



'iKAa 
CeKBa 

^eK'Aa 
\ CeK'Ba / 



in BLG. 



(13) 



(14) 



and use the identifications (12) and (131 of jz and sublattice degrees of freedom. 
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Equations (10) and (lllj) present the Kondo Hamiltonian Hq + Hj in the basis of the operators c^xj^a- This 



Hamiltonian can be rewritten in terms of the chiral ID electrons introduced as 



de 
2tt 



where s < and s > regions correspond to the incoming and outgoing radial waves (we preserve the same notation 
for fields ip{r) in 2D and ip{s) in ID, since this should not lead to confusion). Leaving only the four channels that 
couple to the impurity, the kinetic energy ( 10 ) takes the form of Eq. Q of the Main Text in this basis. The exchange 
interaction (11) takes the local form of Eq. ^ of the Main Text in terms of ipis) if JP'°'{e,e') can be treated as 
energy- independent. This holds exactly in BLG (since its density of states = |^ is e- independent) and holds in 
MLG in the weak-coupling limit, where the value J^'"^ — J''''^{eF, ef) at the Fermi level is taken. 



SUPPLEMENTARY MATERIAL: DETAILS OF THE DERIVATION OF THE KONDO MODELS 

Here, we provide some details of the group-theoretical derivation of the Kondo models summarized in Tab. |ll] of 
the Main Text. 

Impurity orbital singlet class. In the orbital singlet class, Ri = Ai, since the square of any ID IR is a unity IR, 
Ri X rI — Ai/A'^/Aig, the transformation properties of JP'"^ ^ i?o x i?| are determined by those of the conduction 
electrons. 

If in the decomposition of Re (Tab. |l| all IRs are of different type, then only the squares of each IR contain a 
unity IR and produce an invariant, which is a unity matrix in the subspace of the IR. As a result, there are no 
processes transferring electrons between different channels, neither within one 2D IR nor between different IRs: JP'"^ 
are diagonal matrices and within each 2D IR both channels have the same coupling constant, protected by symmetry. 
This is so for classes (I) and (II) but not for (III), i.e., in all cases we consider, except for the BLG (c, Cs^) case. 
Thus, the orbital-singlet class (I) and (II) Kondo models contain two and one symmetric two-channel contributions, 
respectively, decoupled from other channels, a prerequisite for the realization of the two-channel Kondo effect. 

On the other hand, in class (III), the conduction sea R^ = 2E = + E^^ consists of two 2D IRs E"'^ of the 
same type E (the labels a,/3 are used to distinguish between the two subspaces). As a result, the decomposition 
i?e X i?| = 4(Ai + A2 + E) contains four invariants (4Ai): two correspond to the symmetric coupling within each 2D 
IRs and and two, described by one comp/ea; coupling constant J^"^'^ , correspond to "conversion" processes E°' f^- E^ , 
whereby the conduction electrons are transferred between two IRs. 

Impurity orbital doublet class. In the orbital doublet class, Ri = Ei, for any case of the impurity placement, the 
product Ri X Rj = Ai[Ttj] + A2[Tz] + E[Tx^y] decomposes into two ID and one 2D IRs, whose bases are formed by the 
impurity isospin operators indicated in brackets. Likewise, for any given 2D IR Ee in Re, the bases of the IRs in the 
decomposition EeX Eg — Ai[to] + A2[tz] + E[Tx,y] are formed by the unity and Pauli matrices r^, 7 = 0, x, j/, z acting 
in the space of i?e± conduction electron states. (Here, the MLG (5, Csv) case is provided as an illustrative example, 
all other cases are shown in the fourth column of Tab. |l]). 

Evaluating the contributions to the product ([6]), we find that, for any class (I), (II), or (HI) and any 2D IRs Ee 
and Ei, the product {Ee x Ee) x (£"; x £';) contains three invariants, which give rise to an anisotropic orbital Kondo 
interaction Jq TqTq + J'^'TzTz + J'^'ir^Tx + TyTy) between the impurity isospin and the conduction electrons in two 
Ee channels, in both charge {JP) and spin ( J'^) sectors. 

Further, the product {Eei x Ee2) x {Ei x Ei) with different 2D IRs £"01 and £'02 contains no invariants. Thus, in 
class (I), there are no conversion processes Eei -H> Ee2 and the Kondo model consists of two decoupled two-channel 
contributions. In class (II), the products (^1 x A2)e x {E x E)i and (^1,2 x E)e x {E x E)i contain invariants, 
and hence conversion processes Ai o A2, Ai^2 ^ E between all IRs are present. In class (III), the invariants of 
{E X E)e X {E X E)i also result in conversion processes E" ^ E^ , with the structure of the orbital Kondo interactions 
above. , 



